Create new seurat object
norm_mat <- (2**assays(sce)$limma)-1
ln_mat <- log(norm_mat + 1)
srt <- CreateSeuratObject(counts = ln_mat, min.cells = 0, min.features = 0)
## Warning: Non-unique features (rownames) present in the input matrix, making
## unique
## Warning: Non-unique cell names (colnames) present in the input matrix, making
## unique
srt
## An object of class Seurat
## 16864 features across 12313 samples within 1 assay
## Active assay: RNA (16864 features, 0 variable features)
srt <- FindVariableFeatures(srt, selection.method = "vst", nfeatures = 3000)
# Identify the 20 most highly variable genes
top20 <- head(VariableFeatures(srt), 20)
# plot variable features with and without labels
plot1 <- VariableFeaturePlot(srt)
plot2 <- LabelPoints(plot = plot1, points = top20, repel = TRUE)
## When using repel, set xnudge and ynudge to 0 for optimal results
plot2
## Warning: Transformation introduced infinite values in continuous x-axis
## Warning: Removed 49 rows containing missing values (geom_point).
all.genes <- rownames(srt)
srt <- ScaleData(srt, features = all.genes)
## Centering and scaling data matrix
srt <- RunPCA(srt, features = VariableFeatures(object = srt), npcs = 100)
## PC_ 1
## Positive: Itm2b, Apoe, C1qc, C1qb, Selenop, Ms4a7, C1qa, Timp2, Pltp, Trem2
## Hexa, Hexb, Tmem86a, Rgs10, Arhgap22, Tmcc3, Ctsh, Ctsb, Tnfaip8l2, Frmd4b
## Syngr1, Fcrls, Aif1, Serpinb6a, Mertk, Otulinl, Irf8, Itga9, Cxcl16, Maf
## Negative: Clec4e, Pim1, Fn1, Clec4d, Emilin2, Cd14, Ahnak, S100a6, Ehd1, Rab11fip1
## Cebpb, Txnrd1, S100a4, Gda, F10, Pde4b, Nfkb1, Lilrb4a, Ppp4r2, Thbs1
## Cd44, Malt1, Slc7a11, Slfn2, S100a10, Sod2, Adora2b, Plcb1, Ccl9, Srgn
## PC_ 2
## Positive: Ctsd, Fabp5, Ctsl, Cd63, Spp1, Ftl1, Ctsb, Lhfpl2, Pf4, Gpnmb
## Lgals1, Vat1, Trem2, Emp1, Cd68, Fth1, Arhgap10, Basp1, Gde1, Cstb
## Slc6a8, Syngr1, Adam8, Cd36, Lgals3, Il7r, Cpd, Abca1, Fnip2, Timp2
## Negative: Ly6e, Ifi209, Ms4a4c, Samhd1, Plac8, Mndal, Ifi213, A330040F15Rik, Ifitm3, Ifi203
## Txnip, Stat2, Ifi204, Irf7, Rnf213, Dleu2, Slfn1, Ccr2, Ifi211, Slfn5
## Ifi47, Trim30a, Parp14, Isg15, Phf11b, Ifitm6, Plbd1, Arhgap15, Cd74, Ccnd3
## PC_ 3
## Positive: Cytip, Xylt1, Fyn, Runx3, Napsa, Csf2rb, Plbd1, Lsp1, Msrb1, Gngt2
## Itgal, Gan, Cd52, Mir142hg, Hp, Gcsh, Gsr, Antxr2, Gpr141, Pde3b
## Runx2, Tspan13, Lmnb1, Itgax, Cbfa2t3, Gpr132, Mcemp1, Arl4c, Actg1, Gcnt2
## Negative: Tnf, Cxcl2, Ccl2, Ccl7, Nfkbiz, Il10, Tnfaip3, Cxcl1, Nfkbia, Phlda1
## Gm15726, Egr1, Marcksl1, Mrc1, Sdc4, Gadd45b, Sertad2, Mir155hg, Tnfsf9, Sash1
## Errfi1, Orai2, G530011O06Rik, Arhgap10, Maff, Gm14636, Mthfs, Pf4, Ccl12, Nrp2
## PC_ 4
## Positive: Zswim6, Pid1, Trps1, Camk1d, Cmss1, Zfp608, Ssh2, Nedd9, Dleu2, Pecam1
## Grk5, Gm26917, Foxp1, Ints6, Runx1, Zfp710, Fbxl18, Pde4c, Gm26887, Ccnb1ip1
## Slc9a9, Neat1, Sgms1, Dand5, A530013C23Rik, Pdpk1, Ptprj, Ccnd3, Slc12a6, Tmem164
## Negative: Isg15, Rsad2, Irf7, Bst2, Phf11b, Ifit3, Oasl2, Ifi211, Zbp1, Ifit2
## Phf11d, Slfn5, Ms4a4c, Oas3, Ifi47, Ifi204, Usp18, Ifit1, Lgals3bp, Trim30a
## Ifi209, Mndal, Znfx1, A330040F15Rik, Ifi206, Stat1, Parp14, Irgm1, Isg20, Ifi203
## PC_ 5
## Positive: Nfe2l2, H2-Aa, Tmem176b, H2-Ab1, Rel, Il1b, Cd74, Runx3, Ccrl2, Pde4b
## Crem, H2-Eb1, Trf, Tmem176a, Jmy, Cd83, Nr4a3, Lpcat2, Abca1, H2-DMb1
## Tspan13, Fcgr2b, H2-DMa, Kdm6b, Arl5c, Apoe, Itga9, Gpr132, Cxcl16, Mpp7
## Negative: Capg, Anxa2, S100a10, S100a6, S100a4, Lgals1, Anxa1, Crip1, Lpl, Vim
## Lgals3, Pkm, Tagln2, Calr, Rnh1, Pdia6, Tmsb10, Msr1, Adam8, Ccl9
## S100a11, Qpct, Cdk14, Cd180, Cd36, Cd68, Psd3, Hsp90b1, Vat1, Slfn4
ElbowPlot(srt, ndims = 100) +
geom_vline(xintercept = 25, color = "red", alpha = 0.5) +
theme(text = element_text(size = 10),
axis.text.x = element_text(size = 10),
axis.text.y = element_text(size = 10))
srt$sample <- sce$sample
srt$lane <- sce$lane
srt$sorting_day <- sce$sorting_day
srt$cell_type <- sce$cell_type
DimPlot(srt, reduction = "pca", group.by = "cell_type")
srt <- RunUMAP(srt, n.components = 10, features = VariableFeatures(srt))
## Warning: The default method for RunUMAP has changed from calling Python UMAP via reticulate to the R-native UWOT using the cosine metric
## To use Python UMAP via reticulate, set umap.method to 'umap-learn' and metric to 'correlation'
## This message will be shown once per session
## 17:06:23 UMAP embedding parameters a = 0.9922 b = 1.112
## 17:06:23 Read 12313 rows and found 3000 numeric columns
## 17:06:23 Using Annoy for neighbor search, n_neighbors = 30
## 17:06:23 Building Annoy index with metric = cosine, n_trees = 50
## 0% 10 20 30 40 50 60 70 80 90 100%
## [----|----|----|----|----|----|----|----|----|----|
## **************************************************|
## 17:06:33 Writing NN index file to temp file /tmp/RtmpaxeuJs/filed94337e35537
## 17:06:33 Searching Annoy index using 1 thread, search_k = 3000
## 17:08:59 Annoy recall = 100%
## 17:09:00 Commencing smooth kNN distance calibration using 1 thread with target n_neighbors = 30
## 17:09:01 Initializing from normalized Laplacian + noise (using irlba)
## 17:09:02 Commencing optimization for 200 epochs, with 667944 positive edges
## 17:09:12 Optimization finished
srt <- RunTSNE(srt,
dim.embed = 3,
dims = 1:25)
saveRDS(srt, "/mnt/nmorais-nfs/marta/pB_joana/pC_data/srt-mono-macs-young-without-residents-bec.rds")
table(srt$cell_type,srt$sample)
##
## yg0 yg1 yg3 yg5
## CD11a Myeloid 89 7 9 6
## Infiltrating Type 1 293 1371 314 46
## Infiltrating Type 2 26 1865 46 12
## Macrophages intermediate 15 25 2784 2131
## Macrophages Repair 15 12 1254 1922
table(srt$sorting_day,srt$sample)
##
## yg0 yg1 yg3 yg5
## day1 0 3280 4407 0
## day2 438 0 0 4117
table(srt$lane,srt$sample)
##
## yg0 yg1 yg3 yg5
## lane1 0 3280 0 4117
## lane2 438 0 4407 0
table(srt$cell_type,srt$sorting_day)
##
## day1 day2
## CD11a Myeloid 16 95
## Infiltrating Type 1 1685 339
## Infiltrating Type 2 1911 38
## Macrophages intermediate 2809 2146
## Macrophages Repair 1266 1937
table(srt$cell_type,srt$lane)
##
## lane1 lane2
## CD11a Myeloid 13 98
## Infiltrating Type 1 1417 607
## Infiltrating Type 2 1877 72
## Macrophages intermediate 2156 2799
## Macrophages Repair 1934 1269
colors <- c("darkorange4","red", "yellow", "green4", "blue" )
DimPlot(srt, reduction = "tsne", group.by = "cell_type", cols = colors, pt.size = 0.5)
DimPlot(srt, reduction = "tsne", group.by = "sample", pt.size = 0.5)
DimPlot(srt, reduction = "tsne", group.by = "lane", cols = c("orange", "blue"), pt.size = 0.5)
DimPlot(srt, reduction = "tsne", group.by = "sorting_day", cols = c("orange", "blue"), pt.size = 0.5)
DimPlot(srt[,srt$lane == "lane1"], reduction = "tsne", group.by = "lane", cols = c("orange", "blue"), pt.size = 0.5)
DimPlot(srt[,srt$lane == "lane2"], reduction = "tsne", group.by = "lane", cols = c("blue"), pt.size = 0.5)
DimPlot(srt[,srt$sorting_day == "day1"], reduction = "tsne", group.by = "sorting_day", cols = c("orange", "blue"), pt.size = 0.5)
DimPlot(srt[,srt$sorting_day == "day2"], reduction = "tsne", group.by = "sorting_day", cols = c( "blue"), pt.size = 0.5)
DimPlot(srt, reduction = "umap", group.by = "cell_type", cols = colors, pt.size = 0.5)
DimPlot(srt, reduction = "umap", group.by = "sample", pt.size = 0.5)
DimPlot(srt, reduction = "umap", group.by = "lane", cols = c("orange", "blue"), pt.size = 0.5)
DimPlot(srt, reduction = "umap", group.by = "sorting_day", cols = c("orange", "blue"), pt.size = 0.5)
DimPlot(srt[,srt$lane == "lane1"], reduction = "umap", group.by = "lane", cols = c("orange", "blue"), pt.size = 0.5)
DimPlot(srt[,srt$lane == "lane2"], reduction = "umap", group.by = "lane", cols = c("blue"), pt.size = 0.5)
DimPlot(srt[,srt$sorting_day == "day1"], reduction = "umap", group.by = "sorting_day", cols = c("orange", "blue"), pt.size = 0.5)
DimPlot(srt[,srt$sorting_day == "day2"], reduction = "umap", group.by = "sorting_day", cols = c( "blue"), pt.size = 0.5)
DimPlot(srt[,srt$sample == "yg0"], reduction = "umap", group.by = "sample", cols = c( "tomato"), pt.size = 0.5)
DimPlot(srt[,srt$sample == "yg1"], reduction = "umap", group.by = "sample", cols = c( "green3"), pt.size = 0.5)
DimPlot(srt[,srt$sample == "yg3"], reduction = "umap", group.by = "sample", cols = c( "dodgerblue"), pt.size = 0.5)
DimPlot(srt[,srt$sample == "yg5"], reduction = "umap", group.by = "sample", cols = c( "purple"), pt.size = 0.5)
sessionInfo()
## R version 4.1.2 (2021-11-01)
## Platform: x86_64-pc-linux-gnu (64-bit)
## Running under: Ubuntu 18.04.6 LTS
##
## Matrix products: default
## BLAS: /usr/lib/x86_64-linux-gnu/blas/libblas.so.3.7.1
## LAPACK: /usr/lib/x86_64-linux-gnu/lapack/liblapack.so.3.7.1
##
## locale:
## [1] LC_CTYPE=en_US.UTF-8 LC_NUMERIC=C
## [3] LC_TIME=en_US.UTF-8 LC_COLLATE=en_US.UTF-8
## [5] LC_MONETARY=en_US.UTF-8 LC_MESSAGES=en_US.UTF-8
## [7] LC_PAPER=en_US.UTF-8 LC_NAME=C
## [9] LC_ADDRESS=C LC_TELEPHONE=C
## [11] LC_MEASUREMENT=en_US.UTF-8 LC_IDENTIFICATION=C
##
## attached base packages:
## [1] stats4 stats graphics grDevices utils datasets methods
## [8] base
##
## other attached packages:
## [1] SingleCellExperiment_1.16.0 SummarizedExperiment_1.24.0
## [3] Biobase_2.54.0 GenomicRanges_1.46.1
## [5] GenomeInfoDb_1.30.1 IRanges_2.28.0
## [7] S4Vectors_0.32.4 BiocGenerics_0.40.0
## [9] MatrixGenerics_1.6.0 matrixStats_0.62.0
## [11] ggplot2_3.3.6 SeuratObject_4.1.3
## [13] Seurat_4.1.0
##
## loaded via a namespace (and not attached):
## [1] plyr_1.8.7 igraph_1.3.5 lazyeval_0.2.2
## [4] sp_1.5-0 splines_4.1.2 listenv_0.8.0
## [7] scattermore_0.8 digest_0.6.30 htmltools_0.5.3
## [10] fansi_1.0.3 magrittr_2.0.3 tensor_1.5
## [13] cluster_2.1.4 ROCR_1.0-11 globals_0.16.1
## [16] spatstat.sparse_3.0-0 colorspace_2.0-3 ggrepel_0.9.1
## [19] xfun_0.34 dplyr_1.0.10 RCurl_1.98-1.9
## [22] jsonlite_1.8.3 progressr_0.11.0 spatstat.data_3.0-0
## [25] survival_3.4-0 zoo_1.8-11 glue_1.6.2
## [28] polyclip_1.10-4 gtable_0.3.1 zlibbioc_1.40.0
## [31] XVector_0.34.0 leiden_0.4.3 DelayedArray_0.20.0
## [34] future.apply_1.9.1 abind_1.4-5 scales_1.2.1
## [37] DBI_1.1.3 spatstat.random_3.0-1 miniUI_0.1.1.1
## [40] Rcpp_1.0.9 viridisLite_0.4.1 xtable_1.8-4
## [43] reticulate_1.26 spatstat.core_2.4-4 htmlwidgets_1.5.4
## [46] httr_1.4.4 RColorBrewer_1.1-3 ellipsis_0.3.2
## [49] ica_1.0-3 farver_2.1.1 pkgconfig_2.0.3
## [52] sass_0.4.2 uwot_0.1.14 deldir_1.0-6
## [55] utf8_1.2.2 labeling_0.4.2 tidyselect_1.2.0
## [58] rlang_1.0.6 reshape2_1.4.4 later_1.3.0
## [61] munsell_0.5.0 tools_4.1.2 cachem_1.0.6
## [64] cli_3.4.1 generics_0.1.3 ggridges_0.5.4
## [67] evaluate_0.17 stringr_1.4.1 fastmap_1.1.0
## [70] yaml_2.3.6 goftest_1.2-3 knitr_1.40
## [73] fitdistrplus_1.1-8 purrr_0.3.5 RANN_2.6.1
## [76] pbapply_1.5-0 future_1.28.0 nlme_3.1-160
## [79] mime_0.12 compiler_4.1.2 rstudioapi_0.14
## [82] plotly_4.10.0 png_0.1-7 spatstat.utils_3.0-1
## [85] tibble_3.1.8 bslib_0.4.0 stringi_1.7.8
## [88] highr_0.9 lattice_0.20-45 Matrix_1.5-1
## [91] vctrs_0.5.0 pillar_1.8.1 lifecycle_1.0.3
## [94] spatstat.geom_3.0-3 lmtest_0.9-40 jquerylib_0.1.4
## [97] RcppAnnoy_0.0.19 data.table_1.14.4 cowplot_1.1.1
## [100] bitops_1.0-7 irlba_2.3.5.1 httpuv_1.6.6
## [103] patchwork_1.1.2 R6_2.5.1 promises_1.2.0.1
## [106] KernSmooth_2.23-20 gridExtra_2.3 parallelly_1.32.1
## [109] codetools_0.2-18 MASS_7.3-58.1 assertthat_0.2.1
## [112] withr_2.5.0 sctransform_0.3.5 GenomeInfoDbData_1.2.7
## [115] mgcv_1.8-41 parallel_4.1.2 grid_4.1.2
## [118] rpart_4.1.19 tidyr_1.2.1 rmarkdown_2.17
## [121] Rtsne_0.16 shiny_1.7.2